#!/usr/bin/env python
import sys
import itk	
import imageTools.ITKUtils.io as io
import imageTools.ITKUtils.ifilter as ifilter
import numpy as np
import os
import pylab as pb

# get the threshold and seed point
def getThreshSeed(filtered_file, original_file):
    img = io.readImage(filtered_file)
    array = itk.PyBuffer.IF3.GetArrayFromImage(img)
    #pb.hist(array.flat, 100)
    #pb.show()
    mask = np.where( array > 0, 1,0).astype(np.uint8)
    img = io.readImage(original_file) 
    img_array_1 = itk.PyBuffer.IF3.GetArrayFromImage(img)
    img_array = mask*img_array_1
    flat_array = img_array.flat
    index = np.where((flat_array > 100) & (flat_array<600))
    index = index[0]
    flat_array = flat_array[index]
    #pb.hist(flat_array,100)
    #pb.show()
    #pb.hist(img_array_1.flat, 1000)
    #pb.show()
    mean = np.mean(flat_array)
    std = np.std(flat_array)
    value = mean - 0.25*std 
    print mean
    seed = np.where( (array > mean) & ((img_array_1 > value  )& (img_array_1 < 400)))
    seed = np.transpose(seed)
    print seed
    return [value, seed] 

# count the number of voxels included in the segmentation
def count(filename):
    img = io.readImage(filename)
    array = itk.PyBuffer.IUC3.GetArrayFromImage(img)
    count = np.sum(array)
    return count


def main():
   filtered_file = sys.argv[1]
   original_file = sys.argv[2]
   thresh, seed =  getThreshSeed(filtered_file, original_file)
   thresh = int(thresh)
   print seed[0]
   tmp = os.path.splitext(original_file)
   outputfile = str(tmp[0]) + "_seg" + ".mha"
   print outputfile
   xx = int(len(seed)*0.5)
   for j in range(5): 
      print j
      counter = 0 
      for i in range(min(len(seed), 3)): 
         z, y, x = seed[( xx + 13*i + 2*j)%len(seed)]
         cmd = "connectedThreshSeg.py " + "--img "+ str(original_file) + "  --seed  " + str (x) + " " + str(y) + " "+ str(z) + " --thresh "
         cmd =  cmd +  str(value) + " " + str(600)
         print cmd
         os.system(cmd)
         counter =  counter + count(outputfile)
         if counter > 1000:
            j = 10 
            break
      if counter > 100:
            break
      else:
            value = value - 50
 
          
   print counter   
   
   print "done"
if __name__ == '__main__':
    main()


